#include <vector>
#include <cassert>
#include <sstream>
#include <iostream>
#include <fstream>
#include <iomanip>
#include <math.h>
#include <list>
#include <iomanip>
#define M_PI 3.141592653589
using namespace std;
const double EARTH_RADIUS = 6371.004;

extern "C" {
	void get_set(char* path_1, char* path_2, int flag, char* path_3);
	int match(char* path_1, char* path_2);
	double get_out_dis(char* point, char* gis);
}

double rad(double d)
{
	return d * M_PI / 180.0;
}

//计算两个经纬度点之间的实际距离
double GetDistance(vector<double> lg_lt1, vector<double> lg_lt2)
{
	double radLat1 = rad(lg_lt1[1]);
	double radLat2 = rad(lg_lt2[1]);
	double a = radLat1 - radLat2;
	double b = rad(lg_lt1[0]) - rad(lg_lt2[0]);
	double s = 2 * asin(sqrt(pow(sin(a / 2), 2) + cos(radLat1) * cos(radLat2) * pow(sin(b / 2), 2)));
	s = s * EARTH_RADIUS;
	return s * 1000;
}

//转换数组对象格式，获得经纬度对
void get_lng_lat_set(double* m, int len, vector<vector<double>>& lng_lat)
{
	vector<double> temp_lng_lat;
	for (int i = 0; i < len; ++i)
	{
		double lg_lt = *(m + i);
		temp_lng_lat.push_back(lg_lt);
		if ((i % 2 == 1) && (!temp_lng_lat.empty()))
		{
			lng_lat.push_back(temp_lng_lat);
			temp_lng_lat.clear();
		}
	}
}

//获取长路段中距离短路段中某个点最近的点
void get_nearest_point(vector<vector<double>> long_lng_lat, vector<double> lng_lat, size_t& nearest_point)
{
	double dis = 100, temp_dis;
	for (int i = 0; i < long_lng_lat.size(); i++) {
		temp_dis = GetDistance(long_lng_lat[i], lng_lat);
		if (temp_dis < dis)
		{
			dis = temp_dis;
			nearest_point = i;
		}
	}
}

//获取点到直线的距离,此处把任意两个经纬度点确定的线看做直线
double get_dis(vector<double>lng_lat_1, vector<double>lng_lat_2, vector<double>lng_lat_3)
{
	double dis;
	double x;
	double y;
	vector<double> foot;//定义点到直线的垂足
	if (lng_lat_2[0] == lng_lat_1[0]) {//对应斜率为正无穷的情况
		x = lng_lat_1[0];
		y = lng_lat_3[1];
		foot.push_back(x);
		foot.push_back(y);
		dis = GetDistance(foot, lng_lat_3);
	}
	else {
		double k = (lng_lat_2[1] - lng_lat_1[1]) / (lng_lat_2[0] - lng_lat_1[0]);
		double b = (lng_lat_1[1] * lng_lat_2[0] - lng_lat_1[0] * lng_lat_2[1]) / (lng_lat_2[0] - lng_lat_1[0]);
		x = (lng_lat_3[0] + k * lng_lat_3[1] - k * b) / (k * k + 1);
		y = (k * k * lng_lat_3[1] + k * lng_lat_3[0] + b) / (k * k + 1);
		foot.push_back(x);
		foot.push_back(y);
		dis = GetDistance(foot, lng_lat_3);
	}
	return dis;
}

//获取垂足的经纬度
vector<double> get_temp_foot(vector<double>lng_lat_1, vector<double>lng_lat_2, vector<double>lng_lat_3)
{
	double x;
	double y;
	vector<double> foot;//定义点到直线的垂足
	if (lng_lat_2[0] == lng_lat_1[0]) {//对应斜率为正无穷的情况
		x = lng_lat_1[0];
		y = lng_lat_3[1];
		foot.push_back(x);
		foot.push_back(y);
	}
	else {
		double k = (lng_lat_2[1] - lng_lat_1[1]) / (lng_lat_2[0] - lng_lat_1[0]);
		double b = (lng_lat_1[1] * lng_lat_2[0] - lng_lat_1[0] * lng_lat_2[1]) / (lng_lat_2[0] - lng_lat_1[0]);
		x = (lng_lat_3[0] + k * lng_lat_3[1] - k * b) / (k * k + 1);
		y = (k * k * lng_lat_3[1] + k * lng_lat_3[0] + b) / (k * k + 1);
		foot.push_back(x);
		foot.push_back(y);
	}
	return foot;
}
//从长路段中任选两个点确定一条直线，计算任意一点到该直线的距离
double get_point_dis(vector<vector<double>>long_lng_lat, vector<double>lng_lat, size_t ext_point)
{
	double point_dis = 0;
	if (ext_point == 0) {
		point_dis = get_dis(long_lng_lat[ext_point + 1], long_lng_lat[ext_point], lng_lat);
	}
	else if (ext_point == (long_lng_lat.size() - 1)) {
		point_dis = get_dis(long_lng_lat[ext_point - 1], long_lng_lat[ext_point], lng_lat);
	}
	else if (long_lng_lat[ext_point][0] >= lng_lat[0])
	{
		if (long_lng_lat[ext_point - 1][0] < lng_lat[0])
		{
			point_dis = get_dis(long_lng_lat[ext_point - 1], long_lng_lat[ext_point], lng_lat);
		}
		else if (long_lng_lat[ext_point + 1][0] < lng_lat[0])
		{
			point_dis = get_dis(long_lng_lat[ext_point + 1], long_lng_lat[ext_point], lng_lat);
		}
	}
	else if (long_lng_lat[ext_point][0] < lng_lat[0])
	{
		if (long_lng_lat[ext_point - 1][0] > lng_lat[0])
		{
			point_dis = get_dis(long_lng_lat[ext_point], long_lng_lat[ext_point - 1], lng_lat);
		}
		else if (long_lng_lat[ext_point + 1][0] > lng_lat[0])
		{

			point_dis = get_dis(long_lng_lat[ext_point], long_lng_lat[ext_point + 1], lng_lat);
		}
	}
	return point_dis;
}

//获取垂足（求交集并集用）
vector<double> get_foot(vector<vector<double>>long_lng_lat, vector<double>lng_lat, size_t ext_point)
{
	vector<double> point_dis;
	if (ext_point == 0) {
		point_dis = get_temp_foot(long_lng_lat[ext_point + 1], long_lng_lat[ext_point], lng_lat);
	}
	else if (ext_point == (long_lng_lat.size() - 1)) {
		point_dis = get_temp_foot(long_lng_lat[ext_point - 1], long_lng_lat[ext_point], lng_lat);
	}
	else if (long_lng_lat[ext_point][0] >= lng_lat[0])
	{
		if (long_lng_lat[ext_point - 1][0] < lng_lat[0])
		{
			point_dis = get_temp_foot(long_lng_lat[ext_point - 1], long_lng_lat[ext_point], lng_lat);
		}
		else if (long_lng_lat[ext_point + 1][0] < lng_lat[0])
		{
			point_dis = get_temp_foot(long_lng_lat[ext_point + 1], long_lng_lat[ext_point], lng_lat);
		}
	}
	else if (long_lng_lat[ext_point][0] < lng_lat[0])
	{
		if (long_lng_lat[ext_point - 1][0] > lng_lat[0])
		{
			point_dis = get_temp_foot(long_lng_lat[ext_point], long_lng_lat[ext_point - 1], lng_lat);
		}
		else if (long_lng_lat[ext_point + 1][0] > lng_lat[0])
		{

			point_dis = get_temp_foot(long_lng_lat[ext_point], long_lng_lat[ext_point + 1], lng_lat);
		}
	}
	return point_dis;
}

vector<string> stringToken(const string sToBeToken, const string sSeperator)
{
	vector<string> vToken;
	string sCopy = sToBeToken;
	size_t iPosEnd = 0;
	while (true)
	{
		iPosEnd = sCopy.find(sSeperator);
		if (iPosEnd == -1)
		{
			vToken.push_back(sCopy);
			break;
		}
		vToken.push_back(sCopy.substr(0, iPosEnd));
		sCopy = sCopy.substr(iPosEnd + 1);
	}
	return vToken;
}


void txt2vector(char* file, vector<vector<double>> & lng_lat_set)
{
	
	ifstream infile;
	infile.open(file);   //将文件流对象与文件连接起来 
	assert(infile.is_open());   //若失败,则输出错误消息,并终止程序运行 
	vector<string> sub_lng;
	string s;
	const string sep = ",";
	while (getline(infile, s))
	{
		vector<double> lng_lat;
		sub_lng = stringToken(s, sep);
		for (int i = 0; i < 2; i++) {
			double num = 0;
			stringstream ss(sub_lng[i]);
			ss >> num;
			lng_lat.push_back(num);
		}
		lng_lat_set.push_back(lng_lat);
		lng_lat.clear();
	}
	infile.close();
}


//获取某一点到路线的最短距离
double get_point_match_dis(vector<vector<double>>long_lng_lat, vector<double> short_lng_lat) {
	size_t nearest_point;
	get_nearest_point(long_lng_lat, short_lng_lat, nearest_point);
	double dis = get_point_dis(long_lng_lat, short_lng_lat, nearest_point);
	return dis;
}

//获取路线的并、交、补
//线路并集
void get_union_set(vector<vector<double>>lng_lat_set_1, vector<vector<double>> lng_lat_set_2, const char* file)
{
	size_t nearest_point;
	vector<double> temp_point;
	vector<size_t> nearest_points(0);
	vector<vector<double>> union_set;
	for (size_t i = 0; i < lng_lat_set_2.size(); i++) {
		get_nearest_point(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		double dis = get_point_dis(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		if (dis < 1.5)
		{
			if (nearest_points.empty()) {
				nearest_points.push_back(nearest_point);
				nearest_points.push_back(1);
			}
			nearest_points.pop_back();
			nearest_points.push_back(nearest_point);
		}
	}
	if (nearest_points.empty() || nearest_points.size() == 1) {
		vector<double> l(2, 0);
		union_set.push_back(l);
	}
	else {
		if ((nearest_points[0] == 0) && (nearest_points[1] < (lng_lat_set_1.size() - 1))) {
			union_set.insert(union_set.end(), lng_lat_set_2.begin(), lng_lat_set_2.end());
			union_set.insert(union_set.end(), lng_lat_set_1.begin() + nearest_points[1] + 1, lng_lat_set_1.end());
		}
		if ((nearest_points[0] > 0) && nearest_points[1] == (lng_lat_set_1.size() - 1)) {
			union_set.insert(union_set.end(), lng_lat_set_1.begin(), lng_lat_set_1.end());
			union_set.insert(union_set.end(), lng_lat_set_2.begin() + (nearest_points[1] - nearest_points[0]) + 1, lng_lat_set_2.end());
		}
		if ((nearest_points[0] > 0) && (nearest_points[1] < (lng_lat_set_1.size() - 1))) {
			union_set.insert(union_set.end(), lng_lat_set_1.begin(), lng_lat_set_1.end());
		}
		if ((nearest_points[0] == 0) && (nearest_points[1] == (lng_lat_set_1.size() - 1))) {
			union_set.insert(union_set.end(), lng_lat_set_2.begin(), lng_lat_set_2.end());
		}
	}
	ofstream write;
	write.open(file);
	for (int j = 0; j < union_set.size(); j++) {
		write << fixed << setprecision(8) << union_set[j][0] << "," << union_set[j][1] << "\n";
	}
	write.close();
}

//线路交集
void get_intersection_set(vector<vector<double>>lng_lat_set_1, vector<vector<double>> lng_lat_set_2, char* file)
{
	size_t nearest_point;
	vector<double> temp_point;
	vector<vector<double>> intersection_set;
	for (size_t i = 0; i < lng_lat_set_2.size(); i++)
	{
		get_nearest_point(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		double dis = get_point_dis(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		if (dis < 2) {
			temp_point = get_foot(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
			if (intersection_set.empty()) intersection_set.push_back(temp_point);
			else intersection_set.push_back(lng_lat_set_1[nearest_point]);
		}
	}
	intersection_set.push_back(temp_point);
	ofstream write;
	write.open(file);
	write << std::fixed;
	for (int j = 0; j < intersection_set.size(); j++) {
		write << setprecision(8) << intersection_set[j][0] << "," << intersection_set[j][1] << "\n";
	}
	write.close();
}

//线路差集
void get_complement_set(vector<vector<double>>lng_lat_set_1, vector<vector<double>> lng_lat_set_2, char* file)
{
	size_t nearest_point;
	vector<double> temp_point;
	vector<size_t> nearest_points;
	vector<vector<double>> temp_out_set;
	vector<vector<vector<double>>> out_set;
	vector<vector<double>> same_set;
	for (size_t i = 0; i < lng_lat_set_2.size(); i++) {
		get_nearest_point(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		double dis = get_point_dis(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
		if (dis < 1.5)
		{
			if (nearest_points.empty()) {
				nearest_points.push_back(nearest_point);
				nearest_points.push_back(1);
			}
			nearest_points.pop_back();
			nearest_points.push_back(nearest_point);
			temp_point = get_foot(lng_lat_set_1, lng_lat_set_2[i], nearest_point);
			same_set.push_back(temp_point);
			lng_lat_set_1[nearest_point] = temp_point;
			lng_lat_set_2[i] = temp_point;
		}
	}
	if (nearest_points.empty()||nearest_points.size()==1) {
		vector<double> l(2, 0);
		temp_out_set.push_back(l);
	}
	else {
		if (nearest_points[0] == 0 && nearest_points[1] == (lng_lat_set_1.size() - 1)) {
			temp_point.push_back(0.0);
			temp_point.push_back(0.0);
			temp_out_set.push_back(temp_point);
		}
		if (nearest_points[0] == 0 && nearest_points[1] < (lng_lat_set_1.size() - 1)) {
			temp_out_set.insert(temp_out_set.end(), lng_lat_set_1.begin() + nearest_points[1], lng_lat_set_1.end());
		}
		if (nearest_points[0] > 0 && nearest_points[1] == (lng_lat_set_1.size() - 1)) {
			temp_out_set.insert(temp_out_set.end(), lng_lat_set_1.begin(), lng_lat_set_1.begin() + nearest_points[0]);
		}
		if (nearest_points[0] > 0 && nearest_points[1] < (lng_lat_set_1.size() - 1)) {
			temp_out_set.insert(temp_out_set.end(), lng_lat_set_1.begin(), lng_lat_set_1.begin() + nearest_points[0]);
			out_set.push_back(temp_out_set);
			temp_out_set.clear();
			temp_out_set.insert(temp_out_set.end(), lng_lat_set_1.begin() + nearest_points[1] + 1, lng_lat_set_1.end());
		}
	}
	out_set.push_back(temp_out_set);
	ofstream write;
	write.open(file);
	write << std::fixed;
	for (int j = 0; j < out_set.size(); j++) {
		for (int k = 0; k < out_set[j].size(); k++) {
			write << setprecision(8) << out_set[j][k][0] << "," << out_set[j][k][1] << "\n";
		}
	}
	write.close();
}

//获取短路段中所有点到长路段中对应直线的距离的均值以及匹配结果
bool get_match_result(vector<vector<double>>long_lng_lat, vector<vector<double>> short_lng_lat)
{
	vector<size_t>nearest_points;
	vector<double>sum_dis;
	bool flag = false;
	double temp_dis, temp_sum = 0;
	for (int i = 0; i < short_lng_lat.size(); i++) {
		nearest_points.push_back(1);
		get_nearest_point(long_lng_lat, short_lng_lat[i], nearest_points[i]);
		temp_dis = get_point_dis(long_lng_lat, short_lng_lat[i], nearest_points[i]);
		sum_dis.push_back(temp_dis);
	}
	for (int i = 0; i < sum_dis.size(); i++) {
		temp_sum = temp_sum + sum_dis[i];
	}
	double mean_dis = temp_sum / sum_dis.size();
	if ((mean_dis < 2)) {
		flag = true;
	}
	return flag;
}


void get_set(char* file_1, char* file_2, int flag, char* file) {
	vector<vector<double>> lng_lat_set_1;
	vector<vector<double>> lng_lat_set_2;
	txt2vector(file_1, lng_lat_set_1);
	txt2vector(file_2, lng_lat_set_2);
	switch (flag)
	{
	case 0:
		get_intersection_set(lng_lat_set_1, lng_lat_set_2, file); break;
	case 1:
		get_union_set(lng_lat_set_1, lng_lat_set_2, file); break;
	case 2:
		get_complement_set(lng_lat_set_1, lng_lat_set_2, file); break;
	default:
		break;
	}
}

int match(char* file_1, char* file_2) {
	vector<vector<double>> lng_lat_set_1;
	vector<vector<double>> lng_lat_set_2;
	txt2vector(file_1, lng_lat_set_1);
	txt2vector(file_2, lng_lat_set_2);
	bool result = get_match_result(lng_lat_set_1, lng_lat_set_2);
	int out = static_cast<int>(result);
	return out;
}

double get_out_dis(char* point, char* gis) {
	vector<vector<double>> lng_lat_set;
	vector<double> lng_lat;
	string l = point;
	vector<string> sub_lng;
	const string sep = ",";
	double num;
	sub_lng = stringToken(l, sep);
	for (int i = 0; i < 2; i++) {
		stringstream ss(sub_lng[i]);
		ss >> num;
		lng_lat.push_back(num);
	}
	txt2vector(gis, lng_lat_set);
	double dis = get_point_match_dis(lng_lat_set, lng_lat);
	return dis;
}
